library(psych)
library(ggpubr)
library(tidyverse)
library(lmerTest)
We will use regions 1 and 2 as control regions since there are related to motor processing. For ROIs we will use regions 5, 6, 7, 8, 9, and 10 since they all overall with crus I and II which have been the main regions implicated in social cogntion. There are also some studies which point to crus II’s involvement in reinforcement learning.
participants <- read.table(file='../participants_good.tsv', sep = '\t', header = TRUE)
participants
all_sub_roi_con <- read.csv('../univariate_roi/all_sub_roi_contrasts.csv')
all_sub_roi_con$subject_id <- all_sub_roi_con$subj
#all_sub_roi_con <- subset(all_sub_roi_con, select = -c(subj))
# Use "adult as group label
all_sub_roi_con$group[all_sub_roi_con$group == 'control'] <- 'college'
all_sub_roi_con
Use all participants as a continuous age cohort
hyp_rois <- c('TPJ_b-social','ATL_b-social','AMG_r-social','PCC_b-social','dmPFC_b-social',
'vmPFC_b-reward',
'region01','region02','region05','region06','region07',
'region08','region09','region10',
'striatum_ventral','striatum_dorsala','striatum_dorsalp')
hyp_rois_names <- c('TPJ','ATL','AMG (r)','PCC','dmPFC',
'vmPFC',
'CB Region 1','CB Region 2','CB Region 5',
'CB Region 6','CB Region 7','CB Region 8','CB Region 9','CB Region 10',
'Ventral Striatum','Anterior \n Dorsal Striatum','Posterior \n Dorsal Striatum')
all_sub_roi_con$roi <- factor(all_sub_roi_con$roi, levels=hyp_rois)
all_tstats_combined <- data.frame(matrix(ncol=6,nrow=0))
colnames(all_tstats_combined) <- c('task','contrast','roi','t_stat','p.adj','estimate')
tasks <- c('mdoors','social')
contrasts <- c('positive_winVlos', 'all_winVlos')
row_n = 1
for (task in tasks) {
# Filter for task data
task_data <- all_sub_roi_con[all_sub_roi_con$task %in% task, ]
for (contrast in contrasts) {
# Filter contrast data
contrast_data <- task_data[task_data$contrast %in% contrast, ]
for (roi in hyp_rois) {
# Filter for ROI data
roi_data <- contrast_data[contrast_data$roi %in% roi, ]
# Calculate t-test
temp_ttest <- t.test(roi_data$contrast_mean, mu=0)
# Fill in dataframe
all_tstats_combined[row_n,] <- c(task,contrast,roi,
temp_ttest$statistic, temp_ttest$p.value,
temp_ttest$estimate)
# Increase counter for row numbers
row_n = row_n+1
}
}
}
Change numbers to numeric
all_tstats_combined <- transform(all_tstats_combined, t_stat = as.numeric(t_stat),
p = as.numeric(p.adj),
estimate = as.numeric(estimate))
adult_data <- all_sub_roi_con[all_sub_roi_con$group %in% 'college', ]
adult_data_pos <- adult_data[adult_data$contrast %in% 'positive_winVlos', ]
adult_data_pos <- adult_data_pos[adult_data_pos$roi %in% hyp_rois, ]
#allg_data$contrast_mean <- allg_data$contrast_mean * -1
adult_data_pos_mdoors <- adult_data_pos[adult_data_pos$task %in% 'mdoors', ]
adult_data_pos_social <- adult_data_pos[adult_data_pos$task %in% 'social', ]
#allg_data$contrast_mean <- allg_data$contrast_mean * -1
adult_data_pos_mdoors.t_tests = adult_data_pos_mdoors %>%
group_by(roi) %>%
summarise(P = t.test(contrast_mean, mu = 0)$p.value,
Sig = ifelse(P < 0.05, "*", ""),
MaxWidth = max(contrast_mean))
stats <- p.adjust(adult_data_pos_mdoors.t_tests$P,method='fdr')
P_adjust <- stats
Sig_adjust = ifelse(P_adjust < 0.005, "***", ifelse(P_adjust < 0.01, "**",
ifelse(P_adjust < 0.05, "*", "")))
adult_data_pos_mdoors.t_tests$Sig_adjust <- Sig_adjust
p <- ggbarplot(adult_data_pos_mdoors, x = 'roi',
y = 'contrast_mean', add = 'mean_se',
position = position_dodge(0.8))
p + geom_text(aes(label = Sig_adjust, y = 0.3), size = 6,
data = adult_data_pos_mdoors.t_tests) +
theme(axis.text.x = element_text(angle =90, vjust = 0.5)) +
xlab("Regions of Interests") + ylab("Mean Response (positive wins > losses)") +
scale_x_discrete(labels=hyp_rois_names)
ggsave('../univariate_roi/adult_pos_winVlos_ttest_bar_mdoors.png')
adole_data <- all_sub_roi_con[all_sub_roi_con$group %in% 'kid', ]
adole_data_pos <- adole_data[adole_data$contrast %in% 'positive_winVlos', ]
adole_data_pos <- adole_data_pos[adole_data_pos$roi %in% hyp_rois, ]
#allg_data$contrast_mean <- allg_data$contrast_mean * -1
adole_data_pos_mdoors <- adole_data_pos[adole_data_pos$task %in% 'mdoors', ]
adole_data_pos_social <- adole_data_pos[adole_data_pos$task %in% 'social', ]
#allg_data$contrast_mean <- allg_data$contrast_mean * -1
adole_data_pos_mdoors.t_tests = adole_data_pos_mdoors %>%
group_by(roi) %>%
summarise(P = t.test(contrast_mean, mu = 0)$p.value,
Sig = ifelse(P < 0.05, "*", ""),
MaxWidth = max(contrast_mean))
stats <- p.adjust(adole_data_pos_mdoors.t_tests$P,method='fdr')
P_adjust <- stats
Sig_adjust = ifelse(P_adjust < 0.005, "***", ifelse(P_adjust < 0.01, "**",
ifelse(P_adjust < 0.05, "*", "")))
adole_data_pos_mdoors.t_tests$Sig_adjust <- Sig_adjust
p <- ggbarplot(adole_data_pos_mdoors, x = 'roi',
y = 'contrast_mean', add = 'mean_se',
position = position_dodge(0.8))
p + geom_text(aes(label = Sig_adjust, y = 0.3), size = 6,
data = adole_data_pos_mdoors.t_tests) +
theme(axis.text.x = element_text(angle =90, vjust = 0.5)) +
xlab("Regions of Interests") + ylab("Mean Response (positive wins > losses)") +
scale_x_discrete(labels=hyp_rois_names)
ggsave('../univariate_roi/adole_pos_winVlos_ttest_bar_mdoors.png')
allg_data_pos <- all_sub_roi_con[all_sub_roi_con$contrast %in% 'positive_winVlos', ]
allg_data_pos <- allg_data_pos[allg_data_pos$roi %in% hyp_rois, ]
#allg_data$contrast_mean <- allg_data$contrast_mean * -1
allg_data_pos_mdoors <- allg_data_pos[allg_data_pos$task %in% 'mdoors', ]
model_allg_pos_mdoors <- aov(contrast_mean ~ group*roi,
data = allg_data_pos_mdoors)
summary(model_allg_pos_mdoors)
## Df Sum Sq Mean Sq F value Pr(>F)
## group 1 0.01 0.00743 0.097 0.755
## roi 16 1.56 0.09771 1.280 0.202
## group:roi 16 1.16 0.07246 0.949 0.512
## Residuals 1003 76.59 0.07636
allg_data_pos_mdoors.t_tests = allg_data_pos_mdoors %>%
group_by(roi) %>%
summarise(P = t.test(contrast_mean ~ group, mu = 0)$p.value,
Sig = ifelse(P < 0.05, "*", ""),
MaxWidth = max(contrast_mean))
stats <- p.adjust(allg_data_pos_mdoors.t_tests$P, method='fdr')
P_adjust <- stats
Sig_adjust = ifelse(P_adjust < 0.005, "***", ifelse(P_adjust < 0.01, "**",
ifelse(P_adjust < 0.05, "*", "")))
allg_data_pos_mdoors.t_tests$Sig_adjust <- Sig_adjust
p <- ggbarplot(allg_data_pos_mdoors, x = 'roi',
y = 'contrast_mean', add = 'mean_se',
color = 'group', position = position_dodge(0.8))
p + geom_text(aes(label = Sig_adjust, y = 0.3), size = 6,
data = allg_data_pos_mdoors.t_tests) +
theme(axis.text.x = element_text(angle =90, vjust = 0.5)) +
xlab("Regions of Interests") + ylab("Mean Response (positive wins > losses)") +
scale_x_discrete(labels=hyp_rois_names)
#p + stat_compare_means(aes(group = group), label = 'p.signif', label.y=0.3) +
# theme(axis.text.x = element_text(angle =45, vjust = 0.5)) +
# scale_x_discrete(labels=hyp_rois_names)
ggsave('../univariate_roi/allp_pos_winVlos_anova_bar_mdoors.png')
p <- ggbarplot(allg_data_pos_mdoors, x = 'roi',
y = 'contrast_mean', add = 'mean_se',
color = 'group', position = position_dodge(0.8))
p + stat_compare_means(aes(group = group), label = 'p.signif', label.y=0.3) +
theme(axis.text.x = element_text(angle =45, vjust = 0.5)) +
scale_x_discrete(labels=hyp_rois_names)
ggsave('../univariate_roi/allp_pos_winVlos_anova_bar_mdoors.png')
adult_data <- all_sub_roi_con[all_sub_roi_con$group %in% 'college', ]
adult_data_all <- adult_data[adult_data$contrast %in% 'all_winVlos', ]
adult_data_all <- adult_data_all[adult_data_all$roi %in% hyp_rois, ]
#allg_data$contrast_mean <- allg_data$contrast_mean * -1
adult_data_all_mdoors <- adult_data_all[adult_data_all$task %in% 'mdoors', ]
adult_data_all_social <- adult_data_all[adult_data_all$task %in% 'social', ]
#allg_data$contrast_mean <- allg_data$contrast_mean * -1
adult_data_all_mdoors.t_tests = adult_data_all_mdoors %>%
group_by(roi) %>%
summarise(P = t.test(contrast_mean, mu = 0)$p.value,
Sig = ifelse(P < 0.05, "*", ""),
MaxWidth = max(contrast_mean))
stats <- p.adjust(adult_data_all_mdoors.t_tests$P,method='fdr')
P_adjust <- stats
Sig_adjust = ifelse(P_adjust < 0.005, "***", ifelse(P_adjust < 0.01, "**",
ifelse(P_adjust < 0.05, "*", "")))
adult_data_all_mdoors.t_tests$Sig_adjust <- Sig_adjust
p <- ggbarplot(adult_data_all_mdoors, x = 'roi',
y = 'contrast_mean', add = 'mean_se',
position = position_dodge(0.8))
p + geom_text(aes(label = Sig_adjust, y = 0.3), size = 6,
data = adult_data_all_mdoors.t_tests) +
theme(axis.text.x = element_text(angle =90, vjust = 0.5)) +
xlab("Regions of Interests") + ylab("Mean Response (all wins > losses)") +
scale_x_discrete(labels=hyp_rois_names)
ggsave('../univariate_roi/adult_all_winVlos_ttest_bar_mdoors.png')
adole_data <- all_sub_roi_con[all_sub_roi_con$group %in% 'kid', ]
adole_data_all <- adole_data[adole_data$contrast %in% 'all_winVlos', ]
adole_data_all <- adole_data_all[adole_data_all$roi %in% hyp_rois, ]
adole_data_all_mdoors <- adole_data_all[adole_data_all$task %in% 'mdoors', ]
adole_data_all_social <- adole_data_all[adole_data_all$task %in% 'social', ]
adole_data_all_mdoors.t_tests = adole_data_all_mdoors %>%
group_by(roi) %>%
summarise(P = t.test(contrast_mean, mu = 0)$p.value,
Sig = ifelse(P < 0.05, "*", ""),
MaxWidth = max(contrast_mean))
stats <- p.adjust(adole_data_all_mdoors.t_tests$P,method='fdr')
P_adjust <- stats
Sig_adjust = ifelse(P_adjust < 0.005, "***", ifelse(P_adjust < 0.01, "**",
ifelse(P_adjust < 0.05, "*", "")))
adole_data_all_mdoors.t_tests$Sig_adjust <- Sig_adjust
p <- ggbarplot(adole_data_all_mdoors, x = 'roi',
y = 'contrast_mean', add = 'mean_se',
position = position_dodge(0.8))
p + geom_text(aes(label = Sig_adjust, y = 0.3), size = 6,
data = adole_data_all_mdoors.t_tests) +
theme(axis.text.x = element_text(angle =90, vjust = 0.5)) +
xlab("Regions of Interests") + ylab("Mean Response (all wins > losses)") +
scale_x_discrete(labels=hyp_rois_names)
ggsave('../univariate_roi/adole_all_winVlos_ttest_bar_mdoors.png')
allg_data_all <- all_sub_roi_con[all_sub_roi_con$contrast %in% 'all_winVlos', ]
allg_data_all <- allg_data_all[allg_data_all$roi %in% hyp_rois, ]
#allg_data$contrast_mean <- allg_data$contrast_mean * -1
allg_data_all_mdoors <- allg_data_all[allg_data_all$task %in% 'mdoors', ]
model_allg_all_mdoors <- aov(contrast_mean ~ group*roi,
data = allg_data_all_mdoors)
summary(model_allg_all_mdoors)
## Df Sum Sq Mean Sq F value Pr(>F)
## group 1 0.07 0.06990 0.815 0.367
## roi 16 0.85 0.05288 0.616 0.873
## group:roi 16 1.28 0.07986 0.931 0.533
## Residuals 1003 86.05 0.08580
allg_data_all_mdoors.t_tests = allg_data_all_mdoors %>%
group_by(roi) %>%
summarise(P = t.test(contrast_mean ~ group, mu = 0)$p.value,
Sig = ifelse(P < 0.05, "*", ""),
MaxWidth = max(contrast_mean))
stats <- p.adjust(allg_data_all_mdoors.t_tests$P, method='fdr')
P_adjust <- stats
Sig_adjust = ifelse(P_adjust < 0.005, "***", ifelse(P_adjust < 0.01, "**",
ifelse(P_adjust < 0.05, "*", "")))
allg_data_all_mdoors.t_tests$Sig_adjust <- Sig_adjust
p <- ggbarplot(allg_data_all_mdoors, x = 'roi',
y = 'contrast_mean', add = 'mean_se',
color = 'group', position = position_dodge(0.8))
p + geom_text(aes(label = Sig_adjust, y = 0.3), size = 6,
data = allg_data_all_mdoors.t_tests) +
theme(axis.text.x = element_text(angle =90, vjust = 0.5)) +
xlab("Regions of Interests") + ylab("Mean Response (positive wins > losses)") +
scale_x_discrete(labels=hyp_rois_names)
#p + stat_compare_means(aes(group = group), label = 'p.signif', label.y=0.3) +
# theme(axis.text.x = element_text(angle =45, vjust = 0.5)) +
# scale_x_discrete(labels=hyp_rois_names)
ggsave('../univariate_roi/allp_all_winVlos_anova_bar_mdoors.png')
allg_data_pos <- all_sub_roi_con[all_sub_roi_con$contrast %in% 'positive_winVlos', ]
allg_data_pos <- allg_data_pos[allg_data_pos$roi %in% hyp_rois, ]
model_allg_pos <- aov(contrast_mean ~ group*task*roi,
data = allg_data_pos)
summary(model_allg_pos)
## Df Sum Sq Mean Sq F value Pr(>F)
## group 1 0.32 0.3166 4.526 0.033509 *
## task 1 0.65 0.6497 9.286 0.002339 **
## roi 16 2.94 0.1837 2.626 0.000433 ***
## group:task 1 0.19 0.1943 2.777 0.095778 .
## group:roi 16 0.76 0.0476 0.681 0.815591
## task:roi 16 0.52 0.0323 0.462 0.964770
## group:task:roi 16 0.91 0.0568 0.812 0.672752
## Residuals 2006 140.35 0.0700
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
p <- ggplot(data = allg_data_pos, aes(x=roi, y=contrast_mean, fill=task, alpha=group))
p + geom_bar(stat = 'summary', fun = 'mean', position=position_dodge()) +
scale_alpha_manual(values = c(0.3, 1)) +
theme_classic() +
theme(axis.text.x = element_text(angle = 90, vjust = 0.5))
ggsave('../univariate_roi/allp_pos_winVlos_anova_bar.png')
allg_data_all <- all_sub_roi_con[all_sub_roi_con$contrast %in% 'all_winVlos', ]
allg_data_all <- allg_data_all[allg_data_all$roi %in% hyp_rois, ]
model_allg_all <- aov(contrast_mean ~ group*task*roi,
data = allg_data_all)
summary(model_allg_all)
## Df Sum Sq Mean Sq F value Pr(>F)
## group 1 0.23 0.2345 3.129 0.07706 .
## task 1 2.57 2.5675 34.257 5.63e-09 ***
## roi 16 2.23 0.1393 1.858 0.02010 *
## group:task 1 0.74 0.7364 9.826 0.00175 **
## group:roi 16 0.76 0.0473 0.631 0.86084
## task:roi 16 0.38 0.0237 0.316 0.99538
## group:task:roi 16 0.70 0.0435 0.580 0.90098
## Residuals 2006 150.35 0.0749
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
p <- ggplot(data = allg_data_all, aes(x=roi, y=contrast_mean, fill=task, alpha=group))
p + geom_bar(stat = 'summary', fun = 'mean', position=position_dodge()) +
scale_alpha_manual(values = c(0.3, 1)) +
theme_classic() +
theme(axis.text.x = element_text(angle = 90, vjust = 0.5))
ggsave('../univariate_roi/allp_all_winVlos_anova_bar.png')
Social Task